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TECHNICAL MEMORANDUM X- 53628 


THE USE OF DSL / 90 FOR THE SIMULATION 
OF STABILIZED PLATFORM DYNAMICS 

SUMMARY 


To get a feel for the power of the digital simulation language DSL/90 
the linearized X-loop of the stabilized platform for the Saturn guidance system 
was simulated by using this language. The computation times and relative 
errors at different step-sizes are measured and compared with each other. 


INTRODUCTION 


The rapid development of computer technology and the creation of new 
engineering oriented languages have obviously established that general purpose 
digital computers are increasingly suitable for simulation of the dynamics of 
large physical systems. At Marshall Space Flight Center (MSFC) . Huntsville, 
Alabama, a combined effort of the Computation Laboratory and Quality and 
Reliability Assurance Laboratory has been undertaken to simulate continuous 
and discrete dynamics of an aerospace vehicle and its ground support equipment 
(GSE) on a large digital computer. This Aerospace Vehicle Simulation ( AVS), 
also called Launch Vehicle Component Level Simulation ( LVCLS) , is primarily 
being developed as a major tool for the checkout of space vehicles; but hopefully 
it will become an essential part of an integrated engineering information system 
which can be used by various laboratories at MSFC for the design, checkout, 
and management of space vehicles . 

In order to facilitate the input of data into the digital computer and have 
an engineer rather than a programmer control the simulation, several simulation 
languages have been developed. One of the most promising simulation languages 
is DSL/ 90 [1-2]. This language allows an engineer to write a program in close 
resemblance to the block diagram which describes mathematically the continuous 
dynamics of the system to be simulated. DSL/90 also allows one to construct 
any functional block as a macro using FORTRAN statements and to unite macro 


statements with regular FORTRAN statements. The user has the choice of six 
standard integration methods to obtain maximum computation speed. However, 
he is free to design his own integration scheme and call it up whenever he needs 
it. 


For a large scale digital simulation such as AVS, it is necessary to make 
use of a simulation language. In order to get a feeling for the computation time 
of DSL/ 90 and to compare it with that of other available languages, it was pro- 
posed to run simulations of a test case using all available integration schemes. 
The linearized X-loop of the stabilized platform for the Saturn V guidance sys- 
tem was chosen as the test case. The block diagram of this control system, 
which consists of four rather complex transfer function blocks, is shown in 
Figure 1. The total system contains seven real and four complex poles and 
eight real and one complex zeros. 


BLOCK 1 BLOCK 2 BLOCK 3 



TRANSFER FUNCTIONS: 


Gi(s) 


0. 007 115 (l + 17 8069 ) (l + 138 . 31 g ) (* + 187. 192 ) S 

( - \ ( - \ ( - \ ( 2 x 0. 377453 S s 2 \ 
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1 + 5499 
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Input signal E = 10 volts step function at t = 0 


FIGURE 1. BLOCK DIAGRAM OF LINEARIZED X-LOOP 
OF STABILIZED PLATFORM 
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POLYNOMIAL TRANSFER FUNCTION GENERATOR 


A preliminary analysis of the block diagram indicated that the standard 
"blocks" present in DSL/ 90 could not handle conveniently all the transfer func- 
tions in the control system. Therefore, it was decided to use the General 
Transfer Function Macro (GTFM-BLOCK) [3] . 

Since this macro accepts as input the coefficients of a transfer function 
in polynomial form, a FORTRAN program was written to automate the conver- 
sion process from the factored form ( zeros and poles) to the polynomial form 
of the transfer functions. 

The FORTRAN program is described in Appendix A. 


DSL/90- SIMULATION PROGRAM 


There are two ways of handling transfer function blocks. Either one 
uses the three standard DSL/90 blocks [ 1] or one applies the General Transfer 
Function Macro (GTFM) Block, a copy of which was obtained from IBM Federal 
Systems Division, Huntsville. Appendix B explains the GTFM-Block in detail. 
Figure 2 shows the breakdown of the various blocks into standard blocks and 
GTFM blocks. 

Initial runs of the complete system indicated that input errors were 
creating an unstable condition. In order to isolate these errors and correct 
them, each block in the system was run individually with a step input using 
both the macro method and the standard DSL/90 blocks. 

Several simulation runs with the following six standard integration 
methods have been made to compare the computation times. 

1. Milne fifth -order predictor-corrector integration routine, 
variable step-size ( MILNE) ; 

2. Runge-Kutta fourth -order, variable step size (RKS) ; 

3. Runge-Kutta fourth-order, fixed step size ( RKSFX) ; 


3 
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FIGURE 2. BREAKDOWN OF BLOCKS OF BLOCK DIAGRAM IN FIGURE 1 INTO 
DSL/9 O-STANDARD AND GTFM BLOCKS 




















4. Simpson Rule, fixed step size (SIMP) ; 

5. Trapezoidal Rule, fixed step size (TRAPZ) ; 

6. Rectangular Rule, fixed step size ( RECT) . 

The relative error criterion was set to 10 -5 , 10 -4 and 10~ 2 for the vari- 
able step-size methods. For the fixed step-size methods the step size was 
adjusted in relation to the minimum step size of Runge-Kutta (variable step 
size), i.e., to At = 2. 5 x 10 -5 and 1 x i0“ 4 sec, respectively. 

There are three phases of operations: 

1. Translation of DSL/90 program into FORTRAN; 

2. Compilation of FORTRAN into machine language; 

3. Execution of machine language program. 

If the integration method is changed no retranslation and recompilation are 
necessary. The binary deck for executing the simulation rim is the same except 
one control card defining the integration method has to be changed. If a rerun 
with different parameters is to be made, the affected parameter cards only have 
to be changed without a retranslation and recompilation being necessary. 

The DSL/90 program "Simulation of Stabilized Platform" is shown in 
Appendix C . 


RESULTS OF SIMULATION 


The input to the stabilized platform is a step function of E = 10 volts. 

Of the six output signals the output signals BETA, FDBACK and OTiMAC are 
chosen to be plotted for the real-time t = 0 to t = 0. 2 sec ( Fig. 3) . A print-out 
is given in Figure 4. Figure 5 shows the characteristics of the individual blocks 
2, 3 and 4, i.e. , their responses to a step-function. 

Various computation times for the six integration methods available in 
DSL/90 have been measured for different step-sizes and error criteria. They 
are tabulated in Tables I through in. All runs were made on the IBM 7094-11 
with 90 kc tape drives located at the Computation Laboratory. 
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TEST OF DSL/90 (ML02) ON X-LOOP OF STABILIZED 
PLATFORM (RKS) 
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FIGURE 3. OUTPUT SIGNALS FDBACK, OTiMAC AND BETA AS 
RESPONSES TO INPUT SIGNAL E, WHICH IS A STEP FUNCTION 
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FIGURE 5. STEP FUNCTION RESPONSES OF BLOCK 2, 
BLOCK 3 AND BLOCK 4 
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TABLE I. TIMES FOR PREPARATION OF SIMULATION 


Read control cards , mount and rewind tapes 

27 

x) 

sec 

Load translator portion of DSL/90 

19 

sec 

Translate DSL/ 90 program to FORTRAN 

12 

sec 

Compile FORTRAN 

17 

sec 

Load simulator portion of DSL/90 

23 

x) 

sec 

Total 

98 

sec 


x) These times must also be taken into account, if no compi- 
lation and translation precede the execution. 

TABLE II. COMPARISON OF ERRORS AND RUN TIMES FOR 
VARIABLE STEP METHODS 


Method / / 

< & / bas< 

/ 

^ f- -y~ 

/ */ f 

Percent Erro 
?d upon final 

C / 

'/f 

r / Pe rcent Error 

values / based upon maximum values 

*/?/////* 

4/ 

if 

Milne 

10“ 5 

7 x 10 -4 

6. 25 x iO -4 

0.01 

<0.005 

<0. 005 

< 2. 52 

58.40 

0. 98 

11.7 


Milne 

10" 4 

l.Gx lO*’ 1 

2x 10“ 3 

0.005 

<0. 005 

<0. 005 

2. 52 

58.40 

59. 02 

6. 5 


Milne 

10~ 2 

3. lxlO- 4 

2x 10" 3 

0. 005 

<0. 005 

<0. 005 

2. 52 

58.42 

59.02 

8. 9 


Runge-Kutta 
( Variable) 

10~ 5 

1 x 10 -4 

1 .03 x 10' 3 

<0. 005 

<0. 005 

<0. 005 

<0.005 

0. 02 

<0. 005 

7. 1 


Runge-Kutta 
( Variable) 

10“ 4 

1. 25x 10" 4 

4. 1 x 10~ 3 

<0. 005 

<0. 005 

<0. 005 

0. 04 

0.42 

0. 09 

2.6 


Runge-Kutta 
( Variable) 

10~ 2 

5 x 10* 4 

4. 3x 10" 3 

<0. 005 

<0. 005 

<0. 005 

0.01 

0. 98 

0. 19 

2. 0 
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TABLE in. COMPARISON OF ERRORS AND RUN TIMES FOR 
FIXED STEP METHODS 


/ / Percent Error / Percent Error / / 

/ / based upon final values / based upon maximum values / / 

/ ' / */ * / * / * / °v * / ' 

Runge-Kutta 
( Fixed) 

2.5X10 -5 







46. 9 


Runge-Kutta 
( Fixed) 

lx icr 4 

0.01 

<0 .005 

<0. 005 

<0. 005 

0.02 

<0. 005 

12. 8 


Simpson 

2.5xiO~ 5 

<0. 005 

<0.005 

<0.005 

<0.005 

<0. 005 

<0. 005 

31.4 


Simpson 

lxl(T 4 

0 . 01 

<0.005 

<0.005 

<0. 005 

0. 04 

0 . 01 

7.5 


Trapezoidal 

2. 5 x 10 -5 

<0.005 

<0. 005 

<0.005 

<0.005 

<0. 005 

<0. 005 

33. 8 


Trapezoidal 

lxlO -4 

0.01 

<0.005 

<0.005 

<0. 005 

0.05 

0.01 

8.6 


Rectangular 

2.5xlO -5 

0.04 

<0.005 

<0.005 

0.06 

0.65 

1.00 

11 . 8 


Rectangular 

lxl(T 4 

0. 14 

<0.005 

<0.005 

0.23 

2.67 

4. 20 

3. 1 



The total time to prepare the simulation, which includes translation of 
the DSL/90 program to FORTRAN and the compilation of the FORTRAN-program, 
amounts to about 100 sec. This time can be cut in half for a DSL/90 program 
which has already been compiled. The results of Runge-Kutta Fixed Step Size 
with At = 2. 5 x 10 -5 sec are considered as the most accurate values and, there- 
fore, are used as reference. 

Table n shows that the Milne method did not lead to a usable solution 
for this system with any of the relative errors specified. This seems to be 
caused by the low order (rectangular) starter used in the present system and 
to the tendency of the Milne method to become unstable if the system contains 
high-frequency terms [5, 6]. However, satisfactory results have been obtained 
with the Milne method by decreasing the maximum step-size (specifying a small 
DELMAX) at the expense of longer execution time. 

If only signal BETA is requested, a larger step-size can be tolerated, 
resulting in decreased execution time (e.g. , 0. 1 percent accuracy in 2.6 sec 
for Runge-Kutta Variable Step) . However, if other signals, such as OT1MAC, 
are requested, the step-size has to be decreased considerably to obtain the 
same accuracy as for BETA. This is caused by the filtering or smoothing 
effect of Block 3 ( Fig. 5) . For this example, Runge-Kutta Variable Step gave 
the most favorable results with respect to both accuracy and computation time. 
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CONCLUSIONS 


The total time to perform a simulation after the mathematical model of 
the physical system has been established, consists of the time for the following 
procedures . 


1. Writing the DSL/90 program and punching it into cards; 

2. Initializing the DSL/ 90 processor, translating the DSL/90 
program into a FORTRAN program and compiling the FORTRAN program; 

3. Executing the compiled program. 

The execution time of about 3 sec is relatively short for the simulation of a 
medium -size control system with an engineering accuracy of about 0. 1 percent. 
The simulation time is sufficient for on-line non-real -time applications , though 
it is by a factor of about 15 too slow for real-time applications. 

The preparation time of about 100 sec for the punched DSL/90 program 
is nearly thirty times the execution time but still shorter than the time required 
to put a patch board on an analog computer and to adjust its potentionmeters . 

If the same system is to be simulated with different parameters, the preparation 
time shrinks to only 50 sec. 

Writing a DSL/90 program requires less programming knowledge than 
writing a FORTRAN program; however, it still requires some programming 
experience and will not always be readily handled by a practical engineer. 
Nevertheless , it is easier to set up a DSL/90 program than an analog program 
since no scaling, no conversion of transfer functions into differential equations , 
and no hardware limitations are encountered. 

As long as the digital simulation is operated in a batch process mode via 
punched cards, the engineer may prefer the analog computer since he has closer 
communications with the machine. Therefore, any digital simulation language 
such as DSL/90 will gain if it is used on-line, preferably with computer graphics. 
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Though the computation times and numerical errors of the different inte- 
gration methods cannot be generalized, they give a feeling for the order of 
magnitude one can expect. In particular , although our results using the Milne 
predictor-corrector were unfavorable, other users have reported very good 
results when using Milne's method for problems with somewhat different charac- 
teristics. 
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APPENDIX A 


POLYNOMIAL TRANSFER FUNCTION GENERATOR* 

FORTRAN IV 

Purpose and Method 

The purpose here is to generate a transfer function polynomial , given 
a factored transfer function. We use a transfer function of the form: 

k( i + s/pj) ( 1 + s/p 2 ) ... ( 1 + s/P m ) 
or 

2 • p 2 

k( 1 + s/pj) ( 1 + s/p 2 ) ... ( 1 + s/P m 2 ) ( 1 + — s + ) . 

p m p m 

This function is input (k, p«, p 2 , . . . , p ) , yielding the equation 

m 

c n + c<s + c 2 s 2 + ... + c s m in the form ( c 0 , c t , . . . , c ) . 
u i * m nx 

The p *s are converted to double precision complex numbers and sub- 
i 

routine GENER is called to generate the polynomial. 

Usage 


Load the binary deck, followed by the data cards ( see Fig. A-i) . 

F0RMAT ( 13) : N0CASE - number of polynomials to be generated. 

FORMAT (13) : m - number of coefficients in this polynomial. 

( If quadratic term appears , replace m by -m. ) 


" JL — — 

Appendix was prepared by R. A. Lewallen, IBM Corporation, Huntsville, 
Alabama, and Alex Zakson, General Electric Company, Huntsville, Alabama. 
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F0RMAT (6E12.6) k - constant multiplier of expression. 

F0RMAT (6E12.6) p.'s - see METHOD for example of these. 

p. 's are input exactly as they appear here. 

Repeat m, k, and p^'s N0CASE times. 

ACCURACY: Output appears in double precision. Up to 8 posi- 

tions may be accurate, depending on accuracy of 
input. 

RESTRICTIONS: m ^ 99 
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$ IBJOB GO 

$ IBFTC MAIN LIST, DECK 
DIMENSION P( 99) 

DOUBLE PRECISION Q( 199) 

1 FORMAT ( 13) 

2 FORMAT (6E12.6) 

3 FORMAT (//14H INPUT DATA IS/8X, 1HC, 13X, 7HP-S ARE/(6E20.8) ) 

4 FORMAT ( 40H COEFFICIENTS APPEAR IN DESCENDING ORDER) 

C NOCASE IS NUMBER OF POLYNOMIALS DESIRED. 

READ (5,1) NOCASE 
DO 500 J = 1 . NOCASE 

C IF NO QUADRATIC TERM APPEARS, LET M = M. 

C IF EQUATION INCLUDES QUADRATIC TERM, LET M = -M. 

READ (5,1) M 
NOTERM = M 
M = IABS( M) 

READ (5,2) C, (P(I) , I = l.M) 

WRITE (6,3) C, (P(I), I = 1 , M) 

WRITE (6.4) 

N = M 

C DOES QUADRATIC TERM APPEAR. 

IF (NOTERM . LT. 0) N = M-2 
IF (N . LE. 0) GO TO 110 
DO 100 I = 1 ,N 
Q( 2*1-1) = -P(I) 

Q(2*I) = 0. 

C *= C/P(I) 

100 CONTINUE 

110 IF (NOTERM . GT. 0) GO TO 190 
P(M-l) = 2. * P(M-l) / P(M) 

P(M) = 1. / P(M)**2 
C = C * P(M) 

C COMPUTE ROOTS OF QUADRATIC. 

DISCR = P(M-l) ** 2 - 4. *P( M) 

IF (DISCR) 25, 75, 75 
25 Q( 2 * M-3) = - P( M-i) / (2. *P(M)) 

Q( 2 *M-1) = Q( 2 *M -3) 

Q( 2 *M -2) = SQRT (ABS (DISCR) ) / (2. *P(M) ) 

Q(2*M) = -Q(2*M-2) 

GO TO 190 

75 ADD = SQRT (DISCR) / ( 2 . * P^l) ) 

START = -P(M-l) / (2. *P(M) ) 

Q( 2 *M - 3) = START + ADD 
Q( 2 *M - 1) = START - ADD 
Q(2*M-2) =0. 

Q( 2 *M) = 0. 

190 Q(2*M+1) = C 

200 CALLGENER (Q,M) 

500 CONTINUE 
STOP 
END 

FIGURE Al. POLYNOMIAL TRANSFER FUNCTION GENERATOR 
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APPENDIX B 


GENERAL TRANSFER FUNCTION MACRO BLOCK 

Purpose 

This appendix contains a description of a macro that extends DSL/90 
block capability to general transfer functions. 

For the general transfer function of the form: 


INPUT 


G(s) 


OUTPUT 


G<s) 


b Y + b S n 1 + . . . + b 2 S + b t 
n+ 1 n 1 1 


a ,S n + a S n * + ...+a ? S + a 1 
n+i n i 1 


where: 


* . 

n+ 1 


* 0 


1 ^ n ^ K . 


A block diagram method described in detail in Reference 4 can be used to reduce 
this transfer function to block notation. The block notation for this transfer 
function is shown in Figure Bl. This block diagram can be incorporated easily 
into a DSL/90 macro. 
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Description 

For user convenience, a simplified macro argument and table are used 
in this macro generation. Also for computer run time minimization, three 
values of K = 5, 10, 30 have been used. These macros are called T5FCN, 
T10FCN, andT30FCN, respectively. 

To utilize this macro efficiently, the user must supply the numerator 
and denominator coefficients in the form of table arrays. The user simply picks 
the proper macro range (K = 5, 10, 30) and orders the coefficient arrays. In 

th 

particular, for the following m numerator polynomial: 

N j S m + N S m_1 + . . . +N,S + N t 
m+1 m 1 1 


the order of the coefficients for the macro would be: 


b i, b 2. b 3 > 


K-m’ K+l-m 


b K-l’ b K’ b K+l 


where: (a) The first (K-m) coefficients are loaded as zero. That is; 


b< = bo = . . 


= b. . = b__ =0 

k - m-i K-m 


(b) 


The last (m+1) coefficients are the specific m degree transfer 
coefficients in increasing order as power of S . That is for this case: 



K-m m+1 


The order of the denominator follows the numerator order exactly. It is to be 
emphasized that this method assumes the degree of the numerator is equal to the 
degree of the denominator. If the numerator degree is less than the denominator 
it may be manipulated by using appropriate numerator coefficients as zero. 
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The macro is specified by calling it as: 



(Nl, Di) 


Yi - T5FCN (Nl, Dl, XI) 

where: ( 1) Nl is a numerator array of ( 5 + 1) elements , 

(2) Dl is a denominator array (5+1) elements, and 
( 3) XI is the single input to this macro block. 

The numerator and denominator arrays are specified in a table as: 

TABLE Nl ( 1 - 6) = b 4 , b 2 , b 3 , b 4 , b 5 , b 6 , ... 

Di ( 1 - 6) = a 4 , a 2 , a 3 , a 4 , a 5 , a 6 . 

The macros for 10 and 30 degree polynomials have the same format as T5FCN. 
The use of this macro is illustrated in Appendix 3. 



FIGURE Bl. BLOCK DIAGRAM FOR GENERAL TRANSFER FUNCTION 
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APPENDIX C 


DSL/ 90- PROGRAM "SIMULATION OF STABILIZED PLATFORM" 

$1EDIT SYSCKl , SRCH 

SIBLDR MAIN 
$ IB LOR CENTRL 
SIBLDR SEQNM 
$IEDIT 
DECK 

, TITLE TEST OF DSL/90 ON X-LOOP OF STABILIZED PLATFORM -- MOD LEVEL 2 

xMACROOUT - T5FCN ( N, D, IN) 

Mi - INTGRL(0. 0, M) 

M2 - INTGRL(0. 0, Mi) 

J M3 - INTGRL(0.0,M2) 

S M4 = INTGRL( 0. 0, M3) 

M5 - INTGRL(0. 0, M4) 

M - ( 1. 0/D( 6) ) - (IN - D(5) - Ml - D{4) * M2 - D(3) - M3 - D(2) - M4 - D(l) * M5) 
OUT - N( 6) - M 4 N (5) - Mi + N(4) - M2 + N(3) * M3 + N(2) - M4 + N(l) * M5 

^ ENDMAC 

✓ STORAG NUMCF1(6), DENCFl(6), NUMCF2(6), DENCF2(6), 

NUMCF4( 6) , DENCF4( 6) 

* - * COMPUTE BLOCK 1 WITH USER-SUPPLIED MACRO * * * 

STEPIN = 10.0 * STEP(0. 0) 

OT1MAC = T5FCN (NUMCF1, DENCF1, STEPIN) 

* * * COMPUTE BLOCK 2 WITH USER-SUPPLIED MACRO * * * 

J INPT2 = OT1MAC - FDBACK 

] OT2MAC - T5FCN (NUMCF2, DENCF2, INPT2) 

- * * COMPUTE BLOCK 3 WITH STANDARD FUNCTIONS * * * 

Cl STAN = INTGRL (0.0, OT2MAC) 

C2STAN = CMPXPL (0.0, 0.0, 0.0, PI, ClSTAN) 

C3STAN - C2STAN * PI * * 2 
BETA = Cl * C3STAN 

* * * COMPUTE BLOCK 4 WITH USER SUPPLIED MACRO * * * 

^ OT4MAC T5FCN (NUMCF4, DENCF4, BETA) 

FDBACK = C2 * OT4MAC 
/ CONST Cl = 48.284, C2 = 572.958 
PARAM PI - 39. 

TABLE NUMCF1(1 - 6) = 0.0, 0. 724115E-02, 0. 497683E-03, 0. 539198E-05, 

0. 157056E-07, 0.0, 

/ DENCFl(i - 6) = 1.0, 0.708728, 0. 453984E-01 , 0. 352518E-03, 

0. 565167E-06, 0. 101711E-08, 

NUMCF2(1 - 6) = 3*0. 0, 1.0, 1.81851E-04, 0.0, 

DENCF2( 1 - 6) = 3*0.0, 1.0, 0. 386583E-03, 0. 323198E-06, 

NUMCF4( 1 - 6) = 0. 941349E-01, 0. 138514E-01 , 0.7188E-03, 
v 0. 177854E-04, 0. 266930E-06, 0. 101712E-08, 

DENCF4( 1 - 6) = 1.0, 0.708728, 0. 453984E-01 , 0. 352518E-03, 

0. 565167E-06, 0. 101711E-08 
' CONTRL FINTIM = 0.2, DELT = 0. 0002 
INTEG MILNE 

RANGE OT1MAC, INPT2, OT2MAC, BETA, OT4MAC, FDBACK, DELT 
J PRINT 0.002, OT1MAC, INPT2, OT2MAC, BETA, OT4MAC, FDBACK, DELT 
S PREPAR 0.0002, BETA, FDBACK 
GRAPH 2., 2., TIME, BETA, FDBACK 

LABEL TEST OF DSL/ 90 (ML02) ON X-LOOP OF STABILIZED PLATFORM (MILNE) 

END 
STOP 
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